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We compare dynamic mean-field and dynamic cavity as methods to describe the stationary states 
of dilute kinetic Ising models. We compute dynamic mean-field theory by expanding in interaction 
strength to third order, and compare to the exact dynamic mean-field theory for fully asymmetric 
networks. We show that in diluted networks the dynamic cavity method generally predicts magneti- 
zations of individual spins better than both first order ( "naive" ) and second order ( "TAP" ) dynamic 
mean field theory. 
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I. INTRODUCTION 



(Classical) statistical mechanical systems in equilibrium are described by the Gibbs measure, which connects the 
propensity of a system to move between two states taken in isolation (the energy differences between these two 
states) to the probability of finding the system in one of the states, when all states are available. This relation is 
normally used to find equilibrium statistics of a system (magnetizations, correlation functions etc.) by sampling a 
dynamics for which the Gibbs measure in a stationary state. A Markov Chain Monte Carlo (MCMC) method, such 
as Glauber dynamics for Ising systems, which we will review briefly below in Sec. [TTJ can work if the sample average 
converges quickly enough to Gibbs measure, and if the quantity to be measured has wide support in phase space. 
Well-known scenarios where this is case are spin systems in the high-temperature phase and when measuring e.g. 
total magnetization. In the low-temperature phase relaxation time to the Gibbs distribution can be very long. On the 
other hand, if the quantity to be measured is e.g. the magnetization of a single spin, then MCMC in a large system 
is slow for the trivial reason that one needs to sweep through all the spins while just being interested in the changes 
in and average over one of them. If the interactions are weak, marginal probability distributions can be computed 
perturbatively in mean-field theory [IH5], which give closed equations for e.g. single-spin magnetizations. For dilute 
systems, where every spin is not connected to most other spins, very powerful message-passing methods have been 
developed by physicists, information theorists and computer scientists over the last two decades to compute marginals 
of Gibbs distributions quickly and accurately |3J [5] . While these cavity equations cannot (in their simplest form) 
deal with the complex phases of random spin systems at low temperature, in suitable scenarios they are much more 
accurate than mean-field theory, and they greatly improve on MCMC for single site magnetization and other local 
properties by substituting a cumbersome sampling by a direct deterministic computation. The cavity method has 
found to have many technological as well as fundamental applications [5HS]- 

The situation is very different for out-of-equilibrium systems, in itself is an extremely broad term covering every- 
thing from macroscopic hydrodynamics (turbulence) |10j and physical and chemical kinetics |11| to interdisciplinary 
applications of statistical physics to neuroscience, population biology, and other fields [HI [12]. We here consider the 
model systems obtained when generalizing the MCMC rules of Ising spin systems (Glauber dynamics) from the equi- 
librium case (symmetric interactions) to a non-equilibrium case (non-symmetric interactions). Such "kinetic Ising" 
models are only conceptual - but tractable - models of real spin systems driven out of equilibrium, and have mainly 
been studied with applications to neuroscience in mind [14H17) . From the mathematical point of view, they are specific 
examples of Markov chains which do not obey detailed balance conditions. In contrast to equilibrium systems, there 
is hence no simple expression for a stationary state akin to a Gibbs measure, but such a state, when it exists, is a 
(complicated) function of all the details of the model. On the other hand, MCMC works as well in such systems as 
for standard equilibrium Ising models, and mean-field theory have been developed up to second order in the inter- 
action strength [18) . This leaves open the case of dilute kinetic Ising models, where in the equilibrium case cavity 
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methods would be preferable. A dynamic cavity method has only very recently developed for majority dynamics [T5] 
and Glauber dynamics 20J and was investigated by us for parallel and sequential update schemes in |21) . 

The dynamic cavity method as outlined in |21j comprises first an ansatz on probability distributions, similar 
to standard cavity, then the Belief Propagation ansatz that cavity distributions factorize, and then also a further 
assumption of Markovianity. As discussed in [22] . the second assumption is exact for fully asymmetric models with 
the parallel update rule. In a more general case, where either the interaction matrix has both a symmetric and an 
asymmetric component, or where the update rule is different, it is however but an approximation. The numerical 
results of 21j, which showed that for some such mixed instances the dynamic cavity is quite accurate, were somewhat 
unexpected. The main motivation for the present paper is therefore to show more systematically in what parameter 
ranges dynamic cavity converges (for these models), where it is accurate compared to MCMC, and to compare its 
predictions to mean field theory. We will show that for dilute kinetic Ising models, dynamic cavity works also for 
the magnetizations of individual spins, and is considerably more accurate than mean-field calculations of the same 
quantities. 

Kinetic Ising models have been studied by other approaches, and we outline them briefly here. When the dis- 
creteness of states is relaxed to a spherical Ising model. Sompolinsky and Zippelius developed a Langevin equation 
formalism |23) . later extended by Crisanti and Sompolinsky to the non-equilibrium case |24j . where several phases 
of these (dynamical) models are outlined. Although pioneering, predicting magnetizations of individual spins is out 
of scope of such methods, as the sphericity approximation has been made. The dynamic replica theory (DRT) has 
been applied to kinetic Ising models [35], which, by the nature of replica theory, however only applies to averages 
over ensembles of models. Sommers developed a path integral formulation for the Glauber dynamics |26j . which was 
at that time only investigated approximately. As an alternative approach to path integral formulation, generating 
functional analysis was developed to study non equilibrium statistical mechanics of disordered systems [27]. It was 
shown by Ncri and Bolle in |20j that at least in some cases, a dynamic cavity analysis explicitly averaged over a random 
ensemble recovers the results of generating functional analysis. Recently, Hertz and Roudi [2SI HI] used generating 
functional analysis to derive mean-field theories, for infinite-range spin glass models. To compare the accuracy of the 
predictions of single-site magnetizations by the dynamic mean-field formula of |28j to the dynamic cavity for dilute 
mixed models was one further motivation for this work. 

The paper is organized as follows. In section [IT] we describe the Glauber dynamics for spin glasses, the model 
which we will study. In Section |III| we discuss two approaches to a dynamic version of the TAP corrections to first- 
order mean-field theory [181 1281 -30 , while in Section IV we derive dynamic cavity equations for diluted networks in 
parallel update. This derivation should be seen as an alternative and (we hope) clearer alternative to |20j and our 
earlier contribution [21j . The main new results of this paper, on the convergence phase of dynamic cavity and on a 
comparison between the predictions of dynamic cavity and mean-field theory to MCMC are presented in Section |V| 
In Section ( VI I we conclude and discuss possible application areas of dynamic cavity. 



II. THE PARALLEL SPIN UPDATE SCHEME IN DILUTE KINETIC ISING MODELS 

The asymmetric dilute Ising model is defined over a set of TV binary variables a = {a±, . . . , <tat}, and an asymmetric 
graph G — (V, E) where V is a set of N vertices, and E is a set of directed edges. To each vertex Vi is associated 
a binary variable Oi. The graphs G are taken from random graph ensembles with bounded average connectivity. 
Following the parametrization of |27j we introduce a connectivity matrix cy, where cy = 1 if there is a link from 
vertex i to vertex j, = otherwise, and matrix elements cy and Cki are independent unless {kl} — {ji}. The 
random graph is then specified by marginal (one-link) distributions 

c c 

P(dj) = jy^l.cy + (1 - ]y0 fo.cy • (1) 

and conditional distributions 

p(cij I Cji) = eS CiiiCjf + (1 - e)p(cij) ■ (2) 

where i,j £ {1, ...,iV} and i < j. In this model the average degree distribution is given by c, and the asymmetry 
is controlled by e G [0,1]. The two extreme values of e give respectively a fully asymmetric network (e = 0), where 
the probabilities of having two directed links between pairs of variables are uncorrelated, and the symmetric network 
(e = 1) where the two links i j and j — >■ i are present or absent together. The parameter set is completed by a 
(real- valued) interaction matrix Jy. Additional assumptions on the «7y, i.e. smallness or that they are random with 
suitable distribution are stated when used. However, for concreteness the reader may in much of this paper think of 
Jij to be independent identically distributed random variables with zero mean and variance - c (Gaussian or binary) 
such that for the fully connected networks (c = N), the interactions scale as the Sherrington-Kirkpatrick model [3"T] . 



3 



The interactions among spins determine the dynamics of system. In the parallel update scheme, which will be 
considered here, at each (discrete) time, all spins are updated according to the Glauber rule 

,. ,x _ J +1 with probability {1 + exp(-2/3 h t (t + l))}" 1 , . 

a ' l[t + L) ~ \ -1 with probability {l+exp(2/3/i i (t + l))}^ 1 W 

where ftj(t) is the effective field acting on spin i at time step t 

M*) = £ JW*-i) + 0i(*) • (4) 

and the parameter /?, analogous to inverse temperature, is a measure of the overall strength of the interactions. The 
notation j e di in ^ and Q indicates all vertices having a direct links to node i and 0j is the (possibly time- 
dependent) external field acting on spin i. In this paper we will adhere to the convention that the interaction indices 
are written in the same order as the temporal order of the interacting spins. Hence we have Jijai(s)aj(s + 1) and 

Jji<Tj(s)(Ti(s + 1). 

The joint probability distribution over all the spin histories p(a(0), . . . , a(t)) has in principle the following simple 
Markov form 

t 

P(a(0) . . . , <?(*)) = J] W[a{s) | h(s)} p(a(0)) (5) 

s = l 

where W is the appropriate transition matrix describing dynamics and updates. If we would have a fully understanding 
of joint probability distribution defined in ([5| we could compute time dependent macroscopic quantities such as 
magnetization and correlations. The evolution of a a single spin is (trivially) defined by summing over the histories 
of all spins except one 

^(0),...,^))= Y, (6) 

which can be further marginalized to the probability of one spin at one time 

PiM*))= Y Pi(ai(0),...,ai(t)) (7) 

CTt(0),...,CT 4 (s-l),CT i (s+l),...,(T j (t) 

and similarly for pairwise joint probability of the histories of two spins Pjj-(<7j(0), . . . <Ji(i), 0j(O), . . . , Cjit)) and 
Pij{ci(s), <Tj-(s')). Consequently, the time evolution of single site magnetization can be obtained from Eq([7]) as 

mi{t) - "iWPifait)) ( 8 ) 



and similarly the correlation functions 



Cij(s,t)= Y °'i( s ) cr i(<)Pij( 'i( S )! a i(*)) ( 9 ) 



Substituting Eq([8]) into dynamics defined in ^ we have 

( 



i(t) = ( tanh(^ Jjiajit - 1) + 6i{t)) j (10) 



where brackets are average with respect to trajectory history. Equation ( 10 1 is exact for the time-dependent magne- 
tization. It is not directly practical, since the marginal over one spin at one time (the magnetization) depends on the 
joint distribution of all the spins influencing it at the previous time, but as we will see in Section |IIIB| it can be used 
as a starting point of a perturbative calculation. 
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III. MEAN-FIELD THEORIES, TAP, AND THE EXPANSION IN SMALL INTERACTIONS 

The mean field theory of spin glass systems started with the Sherrington Kirkpatrick (SK) model [3T]. In this 
model all spins interact with all other spins (infinite-range couplings), which motivates the simplest mean field or 

"naive mean-field" approximation = tanh/3 (y2j JjiTn-j + #t V Shortly afterwards a more accurate mean field 

theory (TAP) was introduced by introducing Onsager reaction for the SK model. This corrects rrij inside the tanh to 
TYij — f3Jijrrii(l — mj) 2 where Jijirii is the field from spin i on spin j and where \jj — — m 2 ) can be interpreted 
as a local susceptibility at spin j pQ. Since in equilibrium Ising J.y = Jji the TAP equilibrium mean field theory is 

hence to, = tanh (PYlj Jji m j + P^i ~ 2 " m i J2j JjiO- — TO j)) ■ ^ s stated in [V, these results can be derived from the 
cavity approach. These can also be derived by observing that in equilibrium a susceptibility is related to a correlation 
by fluctuation-dissipation, and the appropriate correlation was computed by a perturbative argument pQ. For a 
later approach by field-theoretical methods, expanding a functional determinant describing the fluctuations around a 
mean-field stationary point of an action, see e.g. [32] . and references therein. 

In equilibrium Ising systems the naive mean-field and the TAP approximations can further be computed by ex- 
panding the Boltzmann distribution in the interaction strength [33] . To first and second order in interactions this 
result agree with naive mean-field and TAP. 

Recently, a dynamic version of TAP has been derived by Hertz and Roudi [28l [29] by a field-theoretical argument, 
and we show here in Section |III B| below that this also follows from Information Geometry, essentially a systematic 
expansion in interaction strength. For completeness, we will also show that the same dynamic version of TAP follows 
from the "exact mean- field theory" of Mezard and Sakellariou 30 , as already pointed out in (34] . 

Outside equilbrium fluctuation-dissipation does not hold. Conceptually one could therefore say that "dynamic 
TAP" as such is undefined, or, alternatively, that a proper generalization of TAP to a non-equilibrium system should 
be based on fluctuation relations generalizing fluctuation-dissipation theorems [35] (a task we have not attempted to 
carry out). In this paper we take however a more pragmatic approach, and understand "dynamic TAP" to be the 
formulae derived in [551 US] ■ 

Before turning to the technical discussion, let us note that since mean-field and TAP have obvious computational 
advantages, these theories have been applied in much wider settings than in which they have been derived, particularly 
in neuroscience. For a recent review, see [36] and references therein. 



A. Fully asymmetric networks: a reduced theory 



In this section we recall the theory in 



with a view to compute the expansion in small interactions to third 



order. We start by rewriting the exact equation ( 10 1 in the following explicit form: 



E 

(Ti(t),aai{t- 



p((7 9i (t - 1)) O-j(t) 



h0i(t)) 



2 cosh 



0(E 



jedi •'i* 



o-i(t-l) 



»(*)) 



(11) 



where <7Qi is the collection of spins neighboring i with Cji ^ and p(agi(t — 1)) is the corresponding joint probability 
distribution. In a fully asymmetric network, when an interaction coefficient Jji in above is non-zero, then the opposite 
Jij is zero. Each of the spins <jj(t — 1) on the right hand side therefore does not depend directly on spin i on yet 
one time step before, i.e. on <ji(t — 2). Furthermore, the distribution of each of the <jj{t — 1) will in turn depend on 
distributions of other &k(t — 2), but the distribution of these ak(t — 2) do not depend on the q 7 - (t — 1). If there are 
no short paths in the interaction graph between any pairs of spins Uj on the right hand side of (11) except through 
the cavity spin <Ji(t), or if such paths are unimportant, then the spins <jj(t - 



1) will be independent in an asymmetric 
1) acting on ai(t) will be the sum of independent 



network, and the effective field hi(t) 
random variables. 

When the number of interacting spins is large the distribution of hi (t) follows from the central limit theorem 



p(K{t)) 



1 



: exp 



1) and Vi(t) - (h^t)) 2 



2V; (i) 



(12) 



{*?(*)>• 



where (hi(t)) = 9i{t) + V ,. ,.„ •. ,,,„ , 

the thermodynamic limit (N — > oo) is taken at given connectivity c (so that the terms Jjiirijit 
and then c is taken large (so that there are many of them). In general Vi(t) is defined as 



We note that to arrive at this result, first 
1) are independent), 



Vi(t) 



J2 J H J ki [<o-j(* " l)o*(* - 1)) - m j(t - l )m k {t - 1) 



(13) 
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When the additional assumption that the interaction coefficients Jji arc random, independent and evenly distributed 
and small the sum is dominated by the diagonal terms i.e. 



y»(i) = ^4.(l-m»(t-l)) 

jedi 



(14) 



This gives the "exact mean-field" theory of [30] : 



in 



«> = / 



Dx tanh 



(15) 



with the Gaussian measure Dx = ^= e_T - Equation (15) can be iterated starting from some initial condition to get 
all magnetizations at any time, and is exact when the assumptions hold i.e when the network is fully asymmetric, 
when the cavity assumptions hold, when any spin is influenced by a large number of other spins, and when the 
interactions a re r andom, independent, evenly distributed and small. 

To expand (15) in small interactions we introduce Ci(t) = JjiO- — m i(* — I) 2 ) an d ta ^ e an these quantities of 



order e. We have 



tanh \j3 ((hi(t)) + c t {t)x)] = tanh [P{h t {i))\ + xci(t)p(l - tanh 2 [P(hi(t))]) 

x 2 c 2 /? 2 tanh [P{hi{t))\ (1 - tanh 2 [P{hi(t))\) + C(e 3 ) 



(16) 



where every odd term in this expansion will give zero when integrated against a Gaussian measure. Therefore we 
have 



rm(t) = tanh \fi(hi(t))] - p 2 tanh [P(h t (t))} (1 - tanh 2 [^(/i l (t))])c 2 (t) + C(e 4 ) 



(17) 



We would like to write the right hand side of (17) as tanh [/? ((hi(t)) + Ai(t))] + C(e 4 ). A comparison shows that this 
is possible setting Aj(t) = /3tanh [P(hi(t))] cj(tyWe therefore have to fourth order the following functional expression 



rm(t) = tanh [p ((h % (t)) - p tanh [Pih^t))} c?(t))] + 0(e 4 ) 



To first order in e the solution is 



TOj(i) = tanh 



jedi 



0(e 2 



which is the "dynamic naive mean-field" . Inserting this back in ( 18 ) we have "dynamic TAP" of 



rrii (t) = tanh 



p [ J ri m 3 {t - 1) + 6i{t) - P 2 mi (t) J2 Jjii 1 ~ "»i(t ~ I) 2 

jedi I jedi 



C(e 4 ) 



(18) 



(19) 



(20) 



The last term inside the tanh is of order e 2 and a form analogous to the Onsager back-reaction term; there is no third 
order correction in e in this theory. 



B. The Information Geometry viewpoint 



We will now derive the analogues of (19), (20) and a third order term by following the approach of Information 
Geometry [31 [151 [37]. Let <r(0), . . . , a(t) be the time history of a collection of spins. We assume that these spins have 
been generated by a kinetic Ising model with parallel updates and (possibly) time-dependent external fields. The 
joint probability of the history of all the spins is then 

T 

P(a(0), a(T)\0(Q), 0{T), { J,,}) = J] ^exp^t)/^*))^ coahfoft)), K(t) = 0&) + £ J^t - 1). (21) 
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In Information Geometry the space of these model is considered as a manifold with coordinates being the (many) 
parameters 8(0), ... , 8(T), {Jij}- A sub-manifold is the family of independent models 



P ind (?(0),...,a(r)|^ nd (0),...,^ nd (r)) = nn ex PM^))/2cosh(M*)), h i (t) = 8™\t). (22) 

t=l i 



A mean-field approximation is defined as the the independent model with the same magnetizations as the full 
model [37]. For our case it is easily seen that nii(t) = tanh(#j nd (t)) is the variational equation with respect to 
parameter 6>j nd (i) of the Kullback-Leibler divergence Z?_i [p|p md ] = J^P^^a- Therefore, the mean field approxi- 
mation in Information Geometry can also be seen as the independent model with the least Kullback-Leibler divergence 
from the full model [3 EE] [37]. 

Following the approach of |18j we take the interaction parameters (Jij) as small parameters (of order e), and assume 

that the differences A0 4 (*) = 8i(t) - 6»j nd (i) can be expanded in e: 

A8 t (t) = eA^( s ) + e 2 A ( i 2) (s) + ... (23) 
We can then write in analogy with Eq.3.2 of [TH] 

i.s 4 ^ ' j,k 

E tf tndAf ' W + f £ ^l M A"'e J A<"e K + M (24) 

Here 0j stands for the set of all interacting couplings and external fields and J runs over relevant indices. The 
subscript indicates that all derivatives are evaluated at the independent model, and the left-hand side is zero because 
this is the variational equation. In the last term the sum goes over all the parameters labeled J, K and the parameter 
increments are the first order terms A^ (s) and JkV, on third and higher orders mixed terms of A^(s) and A| (s) 
will appear. A calculation presented in Appendix gives the results 



= -E^m^i-l) (25) 

3 

Af\t) = m t (t) Y, - m\(t 1)) (26) 

k 

AfW = -^(l-m 2 fc («-l))J fci Af (t-1) (27) 



Equation ( 23 ) together with the variational equation can be re- written 



tanh- 1 mi (t) = 8 % (t) - eA^(s) - e 2 Af \s) - e 3 Af \s) + 0(e 4 ) (28) 



It is seen that to e this is "dynamic naive mean- field", compare (191, to e this is "dynamic TAP", compare (20), and 
to to e 3 typically there is a non-zero term absent in (20). Such a higher-order difference between the exact mean-field 
theory for the asymmetric model and the field-theoretical approach of 28, 29 was also pointed out in [3U] (page 4, in 
text below Eq. 7) . 



IV. DYNAMIC CAVITY METHOD 



The cavity method for equilibrium systems was introduced in [38, 39 while the dynamic version was studied but 
recently [T9H2T] . In contrast to the equilibrium case, using only the cavity assumption does not in general provide us 
with an efficient algorithm in the dynamic case, but further assumptions are necessary. In this section we derive the 
dynamic cavity method for the kinetic Ising problem, taking a more explicit route than in [20j and [21) . 
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A. Cavity and BP on spin histories 



We consider a number of spins evolving according to a dynamics such as ^ , and we let X denote the whole history 
of spin i, Xi = (<Ti(0), <7,-(l), . . . , cr T (0)). The probability in ^ can then be alternatively be interpreted as a joint 
probability of spin histories, P{Xi 1 X 2 , ■ ■ ■ , X/v), and this probability can be represented by a graph where nodes i 
and j are connected if either or Jji (or both) are non-zero. The corresponding product form is 

P(Xi,X 2 , . . -,X N ) = TJ e ^° 6i W <Ti ( s ' ) T\e^° Ji 3 <?i(s)<? 3 (s+i) TTg-E.iosacoshCfliW+Ej (29) 
which is already normalized. Belief Propagation is expected to work well if this graph is locally tree-like i.e. if all 



loops are long, and can be ignored [5], In (29 1 this is never the case, even if the couplings are fully asymmetric, 
for the simple reason that if spins i and j drive spin k, then they are coupled both by the terms <7i(t)o~k(t + 1) 
and aj(i)<Jk(t + 1), and by the normalization log 2 cosh hk (t + 1). However, these couplings are of a rather peculiar 
type. To proceed we introduce four different marginal probabilities. The first P$ is the marginal probability of spin 
history Xj. The second Pi+di is the marginal probability on the set of histories {Xi{J Xgi}. The third Pgi is the 
marginal on the set of histories Xgi . The fourth and last is PW , a cavity distribution on Xgi . This we take as the 
marginal on Xgi in a revised model where both the spin history Xi as well as the normalization log 2 cosh hi (i) have 
been eliminated. All four probabilities depend on external field parameters which are not necessarily the same. In 
particular, we will express Pgi with one set of external fields as P^> with another set of external fields. By definition 
Pi+Qi — W(Xi\ Xgi)Pgi. The peculiarity of the model is that the (normalized) conditional probability VK(X;| Xgi) is 



already explicitly included in ( 29 1 . We can therefore compare 



W{Xi\Xgi) Al 

TT [ e E s O k (s)a k (s) JJ e J2 B J fc ia*WaK*+l)rr e -E,log2coBh(e*(a)+Ei ^^i(s-l)) j ^ 



to 



P^(Xgi) OC e E 3 e i( s ) ff j( s )+E fce aj,fc^i E s ^fcj CT fc( s ) CT 3( s + 1 )-E s lo g 2cosh ( e 3(s)+E fc ^i •■'fcjCTfc(s-l)) 

X\{X, U^a,} fc #^ \ (safe fc / 

This comparison shows that Pgi with external fields 9j(t) is the same as PW with modified external fields + 
JijO~i{t — 1). Since Pi(Xj) = ^9i)-PaiP^9i) we can therefore write the marginal probability Pj as 

P(X I 0,(0), 9 % (t), ■) = J2 P {l) (Xm(0) |^(0), . . . , 0«(t), •) 

O-8i(0)...CT8i(t-l) 
t 

n^( ff i( s )IM«))ft(M0)) ( 32 ) 

s=l 

where • indicates all the parameters (external fields and interactions) which are the same on the two sides of the 
equation, and 

9f{s) = e j {s) + J i:i a i (s-l) s = 0,...,t (33) 

is the set of external fields that arc modified. 

The next step is to make the Belief Propagation assumption that the spin histories are taken independent in the 
cavity graph: 

p (i) (x 9i \e®(o),...,e$(t),-)=l[ n j ^(x j \9?\a),...,ef\t) i -) (34) 

jedi 
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We here used used /Uj->.j to represent marginal probabilities of neighboring spins, as standard in the BP literature. 
We now consider the subgraph connected to one spin j in the cavity of i. In analogy to above we want to 

compare the marginal on the set of neighbours to spin j in ly , to the cavity distribution on the same set of variables. 
As above the first with one set of external fields is the same as the second with another set of external fields, and we 
therefore find the following recursion equations ("BP update equations"): 

^ l (x ] \ef(V),...,ef(t),- T ^)= II ^ j (xM Ui \^---^ Uj) (t-^- T ^) 



n^ooi^oo) m^(^(o)) 



(35) 



where ■„(<) indicates all the parameters (external fields and interactions) which are the same on the two sides of the 

3 

equation, and 



9 (0. «\ s )=9®(s) + J jk * j (s-l) a = 0,...,t 



is the set of external fields that are modified. Since in fact 6 k l \s) — #fc(s) (spin k is not directly connected to i) we 



(36) 



note that in (36 1 the upper index (i) can be dropped on both sides. The effective field on spin j in Tj is 

h ?( s ) = E JkjVki'-V + Ojia) 

k£dj\i 

and Wj((7j | hj(s)) is the transition probability for the single spin j in the model on 



(37) 



The marginal probability over the history of one spin ( "BP output equation" ) follows from ( 32 1 and ( 34 ) and is 



p i (x i \e i (o) ) ...,e i (t) ) -)= II ^,(^(o),..,^((-i)|f(o) I ..,9f(i-i),.) 

X gi kedi 

t 

11 Wi(o-i(s)|^(s))jii(£ri(0)) (38) 



Equations ( 35 ) and ( 38 1 are the dynamic cavity equations for our system. Both are large sets of equations connecting 



marginal distributions and cavity distributions between two probabilistic models with different parameters. In general 
these equations are (as far as we know) only of conceptual value since on top of connecting different models, the right 
hand side also involves on the order of 2 T ' dl ' operations. In BP such an operation would have to be iterated (for all 
variables) a number of times to reach convergence: as T grows large this becomes unfeasible for the same reason that 
ordinary BP does not work well if the state space of each variable is large. 

We can define (and will later use) marginalizations of the messages down to one time (it is no restriction to take 
this time as the last time): 



ifo(*)i0j o (t))= E 

(T j (0)...f7j (t- 



^i(«T i (0),..., O - J -(t)|df(0) 



(*)) 



(39) 



but in general these quantities do not obey closed equations among themselves. An important exception are fully 
asymmetric networks, since there at most one of Jij and Jji is non-zero. We note that in (35 1 and (38 1 the probability 



distribution of spin i depends on the neighbors di through the effective fields hj(s) and /ij(s), but the messages sent 



from the neighbors to i also depend parametrically on the history of i through the modified external fields 
back-action is absent for the fully asymmetric case where 9^ = Ok independent of spin i. 



7 fc ■ 



This 



B. The projected dynamic BP 



As discussed the marginalization of dynamic BP over one time is not in general a Markov process. However, the 
long time behavior of dynamics (stationary state) is often demanded in many cases. In this section we explain an 
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approximation scheme for computing marginal probabilities for one spin over one time in stationary state, a procedure 
called one-time approximation in 2() and time factorization in [211 122) . 

We start with the dynamic BP for the time histories of messages, i.e. Eq(35l, where on the right hand side the 
time trajectory of messages sent from neighboring spins carry the information from the whole time history of those 
spins. We note that the full dynamics Eq. ^ is in fact Markov. This, and the need to introduce some approximation, 
motivates the time- factorization ansatz which we write for the terms in the right-hand side of Eq(35|: 



z 

Mfe^i (ffc(O), • • ■ ,o-fc(i) |- T o>) = Yl MjUj (o- k {s) \- T U) 



(40) 



where ■ rp <j) indicates the parameters of the model, the same on both sides. Obviously, when inserted into the right- 

1 k 

hand side of (351 such a factorization is not preserved on the left hand side. Since we deal with binary variables we 



rpKJ } 

can introduce time-factorized cavity biases u fc ^-(s), again written for the right-hand side of (35) which are defined by 

(Ws) I - T u 



Mfc- 



2 cosh 



(41) 



A crucial observation is now that when the time-factorization ansatz has been made the cavity biases at different 
external fields are simply related. We will need 



u k^ 3 ( s ) = u kUj( s ) + J jk.o-j{s - 1) s = 0,...,t 



(42) 



which follows from the relation (36). Inserting (41) and (42) into (35) gives 



lij^i{crj{Q),...,o-j{t)\ ■„(*)) 



E 



n n 



/3<x fe (s) «fe^W+Jjk^(s-l) 



x aA< (o),...,<r aAi (t-i) s=0kedj\i 2 cosh 



n(0 



P U kUj( S ) + J jkOj{s - 1) 



(43) 



This equation can be marginalized explicitly over the last time to give 



E 



n 



<r J (t-2),o- eJNi (i-l)fce8j\t 2 cosh 



T («) 



P[u^ {t-l) + J. jk a J (t-2) 



(44) 



The projected dynamic cavity is then to use (44) to compute the terms in a time- factorization of the left-hand 
side of Eq(35). Except for fully asymmetric models (with parallel updates), this approach is not appropriate for 
transients [22 . However, when the external fields 0i are constant in time and when a stationary state has been 
reached, it may be acceptable to also take the messages independent in time. For one and the same set of parameter 
values the fixed-point equations for the time-independent time-factorized cavity biases are then 



1 

w 



E a 3 l0 ' 



E 



aP T.k£dj\i vi* (uj^j +J j k<?' J ) 



i 2 cosh 



2cosh(P hf) 2 cosh ((3u*^) 



,W - 



E J kj°k + 8j 



(45) 



k<£dj\i 



Eq (45) is as ordinary BP solved by iteration, where the right-hand side is computed from u^_^ at iteration time 

t — 1, giving the left hand side at iteration time t. The spin o~j summed over is then conceptually at time t, the 
spins tJfc at time t — 1 and the last spin cr' at time t — 2, all these in the iteration time. 
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Using the iteration time as a proxy for real time we note that in a transient we can compute the time evolution of 
magnetization which would follow from (41 1, (37) 

e /3Efc e 8i\ 3 K^i(t-i)+Jik^(*-2)kk(t-i) 

m (t) = - 

, aiXj( ti),a i{ t-2) Ukeai\j 2 ^sh[p(u k ^ i (t-l) + J ik a i (t-2))] 



tanh 



jedi 



2cosh(/3u i _ >J (< - 2)) 



(46) 



This is not expected to be accurate unless we are already in a stationary state. We use it below in Section |V| as a 
proxy to monitor if the system is in a stationary state. 



V. RESULTS 



In this section we investigate the performance of dynamic cavity method in computing stationary states of diluted 
spin glass in parallel update, and compare to MCMC (Glauber dynamics) and to dynamic mean-field and dynamic 
TAP as defined in Section III The convergence of projected dynamic cavity (dynamic cavity in time-factorized ap- 
proximation) is monitored by comparing magnetization computed from ( 46 ) at successive times for different parameter 
values of the model, and these predictions are then compared to dynamic mean- field and dynamic TAP and MCMC. 



A. Convergence of dynamic BP 



In order to detect where dynamic BP reaches a stationary state we compare single magnetization in two successive 
time step as 

N 

A(t) = 1/NJ2 K(*) - m >tt - !)) 2 (47) 
i=i 

Whenever this deviation vanishes dynamic BP must have converged to a stationary state. Fig. [I] shows the results 
for various connectivity parameters in symmetric and partially symmetric networks. In high temperature we observe 
convergence towards a fixed point whereas in low temperature BP does not reach a stationary state. Roughly speaking, 
dynamic BP stops converging at a value f3 cr (c) which depends on average connectivity. In Fig. ^ the convergence of 




FIG. 1: (Color online) Squared deviation of spin averages between successive update A(t) = 1/iVE^i (™i{t) - m^t - l)) 2 at 
stationary limit for different values of average connectivity c. Mean magnetizations are calculated by projected dynamic cavity 
method i.e Eq. [46] Left panel: partially asymmetric networks with e — 0.5. Right panel : symmetric network e = 1. The 
results are averaged over BP initial conditions (10 experiences). System size is 1000 and external fields are set to zero. 



dynamic BP is plotted to show the effect of asymmetry. In this case it is simply so that for very asymmetric graphs 
BP converges in a very wide region, presumably for arbitrarily large values of (3 if the network grows large enough, 
and, in general, the more asymmetric the network, the better the convergence. 
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<0.4 



H — h £ = 0, C = 2 

H £ = 0.5, c = 2 
Me = 1, c = 2 





<0.5 




FIG. 2: (Color online) Effect of asymmetry (e = 0, 0.5, 1) in squared deviati on of spin averages between successive update 
A(t) = 23j=i ( m »W — m,i(t — l)) 2 , obtained by projected dynamic BP Eq. 46 at stationary limit for average connectivity 



c = 2 (left panel), c = 3 (middle panel), c = 4 (right panel). The results are averaged over BP initial conditions (10 experiences). 
System size is 1000 and external fields are set to zero. 



B. Performance of dynamic BP 



Fig. [3] shows a comparison between dynamic cavity method and dynamic mean field for total magnetization in 
spin glass systems with different asymmetric parameter. The results are obtained in present of small external fields 
6 = 0.001. Dynamic cavity method shows a strong agreement with numerical simulations of type Glauber dynamics 
when it converges to a stationary state. The dynamic mean field method however starts to deviate from numerical 
simulations already in small /3 indicating that it is less accurate compared to the dynamic cavity method. 



0.04 



0.03 - 
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S 

c 0.02 - 
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-i 1 r- 
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#-9 £ = l,c = 3 











FIG. 3: (Color online) Mean square error 8(t) = l/-AT^V =1 (m? 



predicted 



(*) 



empirical 



(t)) 2 of the two approximation methods 



(dynamic mean field and dynamic cavity with respect to the empirical data (Glauber dynamics). Left panel: dynamic mean 
field (Eq. 19 1 for networks with different asymmetric parameter (e = 0, 0.5, 1) and fixed average connectivity c = 3. Right panel 
: the corresponding results obtained by the projected dynamic cavity method Eq. 46 1 For small ft, i.e. high temperature they 
are in agreement with numerical simulations. In low j3 however, dynamic BP outperforms dynamic mean field. 



In order to observe the comparison in more detail, we show also the scatter plot of spin-by-spin magnetization in 
Fig|4] Dynamic cavity method predicts perfectly local magnetizations for fully asymmetric networks and agrees quite 
well with numerical simulations in high temperature for fully symmetric network whereas naive mean field and TAP 
start to deviate already at moderate temperatures. 



VI. CONCLUSION 



Message-passing methods have become an important topic on the border-line between equilibrium statistical physics 
and information theory. In the present paper we have studied an extension of message-passing to non-equilibrium 
Ising spin systems. In contrast to the equilibrium case, the cavity method is not immediately useful to describe the 
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Fully asymmetric e = 





empirical data ( Glauber dynamics) 



1.4 -0.2 0.2 0.4 
empirical data ( Glauber dynamics) 



FIG. 4: (Color online) Scatter plot of local magnetizations for dilute asymmetric networks for four different temperature 
P — 0.4,1,2.5,3.5 and a fixed average connectivity c = 3. Local magnetizations are obtained by dynamic mean field Eq. |19| 
(green), dynamic TAP Eq. 20 (red), and projected dynamic BP Eq. 46 (blue). Left panels show the scatter plots in fully 
asymmetric networks (e = 0) where the projected dynamic BP provides exact results and right panels are scatter plots in fully 
symmetric networks (e = 1). In high temperature, all three methods agree with numerical simulations. In low temperature BP 
starts to outperform naive mean field and TAP. 



dynamics, even if the topology is suitable, because the messages depend on whole spin time histories. The time- 
factorization assumption, as discussed here and in [20H22] . (or some other simplifying assumption) is necessary to 
reduce the complexity, but when so doing one is generally restricted to stationary states. 

We have studied dynamic cavity in the time-factorized assumption for stationary states and outlined its convergence 
region in parameter strength (/?), connectivity (c) and asymmetry (e). By analogy with generally known facts about 
BP it can be argued that when dynamic cavity converges it should typically be a good approximation; the region of 
convergence is therefore a useful proxy for the accuracy. Expanding on first results presented in [21] we show that 
the convergence region in j3 increases with the connectivity. We also find that the convergence region increases with 
asymmetry for several values of connectivity, and that it converges for any interaction strength for fully asymmetric 
networks (as expected) . For networks of moderate size we have directly compared dynamic cavity and dynamic mean- 
field to direct simulation. For several values of asymmetry and connectivity we find that their convergence regions are 
very similar, if not identical, but when both methods converge, then dynamic cavity is considerably more accurate, 
except in the low /3 limit where their performance is about the same. We have hence showed that dynamic cavity 
can be useful new approximation to the dynamics of non-equilibrium spin systems - and any system which can be 
fruitfully modeled by such methods. 

On the analytical side we have discussed the special status of fully asymmetric models, for which the cavity approach 
is in some sense exact. We have also re-derived the "dynamic TAP" equation of Hertz and Roudi [551 HH] using a 
straight-forward approach borrowed from Kappen and Spanjers' treatment of the stationary state |18j clarifying 
that this approach is based on minimizing the distance, in the sense of Information Geometry, to the sub-family of 
independent (but time-changing) models. Whether such a pcrturbative argument can be extended to small deviations 
from e.g. fully asymmetric models remains to be seen. 
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Appendix A: The Information Geometry calculation to second order 

The following calculations are completely parallel to those in Appendix 1 of [IS] and start from 
drrii(t) . 



d0j{a) 
dmj{t) 
dJjk 
d 2 m t (t) 



de 3 {s)de k { s ') 

d 2 mi (t) 



dJ jk d6i(a) 
d 2 m t {t) 



dJjkdJl m 



ind = 5 s,t8ij(l - m%(t)) 
ind = s ikO--m^(t))m j (t-l) 
ind = -2rrii(t)(l - m 2 (t))S ij S lk 5 St tS s ^t 

ind = -2mi(t)(l - m 2 (t))mj(t- l)S ik SuS Stt 

+(1 - m?(t))(l - m 2 (t - l))S ik S kl S Slt -i 
ind = M 1 - »»?(t))(l - m){t - l))5i k mi{t - 2) + (jk) O (Zm) 

-2771^^(5^(1 - m 2 (t))S im (m k (t - l)m m (t - 1) + <W(1 - m^(rj - 1))) 



To first order in e ( 24 ) hence gives 



J2 ~ w'W)Af + J2 M 1 ~ m?(t))ro fc (t - l)J jk = 



which is simply 



4 1) (t) = Af ) (t) + ^J ji m J (t-l) = 



This is the same as "dynamic naive mean field" 

tanh _1 (mi(t)) = 0,(1) + ^ J ki m k (t - 1) + 0(e 2 ) 

k 

The terms arising from second order derivatives and first order increments can be grouped together as 
(1 - m 2 (t)) I - mi (t)(A^) 2 (t) - £(1 - m 2 (t - l))Jii4\t 1) - m*(t) ]T J 2 ki (l - m 2 (t 1)) 



(Al) 
(A2) 
(A3) 

(A4) 
(A5) 
(A6) 

(A7) 

(A8) 



(A9) 



which together with the first order conditions (A7) and the term from the first order derivative and second order 
increment (1 — mf (i))A> 2 '(t) gives 



This is the same as "dynamic TAP" , compare ([20]) above 

tauh _1 (mi(t)) = 9i(t) + ^ Jum k (t - 1) - m;(f) ^ J*rog(t - 1) + G(e 3 



(A10) 



(All) 
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Appendix B: The Information Geometry calculation to third order 

Third order contributions consist partly of terms involving lower than third order derivatives and higher than first 
order increments. The calculation of these use the same elements as above and are 

E^f lindAf( S ) = (l-m^))Af>(*) (Bl) 



j.s.k.s' - 1 



S d jZf(sM JjkA ^ {s) = -2"»i(t)(l-m?(t))X;m J -(t-l)J j <A( 2) (*) 

j,k,l,s jk U ' j 

(1 - m?(t)) - ™ 2 fe (* - l))^A< 2) (t - 1) (B3) 



where two terms can be combined to 



The remainder is 



2m i (t)(l-m 2 i (t))A ( ?\t) ^\t) + £ Jum k {t - 1)^ =0. (B4) 

(1 - m 4 2 (t)) ^Aj 3) (t) + 2(1 - m fe(* - 1)) ■JfciAjj. 2 ^ - 1)^ (lower order terms) (B5) 

To proceed with the terms from third order derivatives and first order increments it is useful to introduce the 
streamlined notation 

mi — nii(t) m[ = m,(t — 1) m'( = rrii(t — 2) etc. (B6) 

and similar for all other quantities. It is also useful to note that though the derivatives act on the complete expres- 
sion involving both probability density P and the tanh they partially obey a chain rule when taken to act on the 
magnetizations alone: 

• a derivative with respect to an external field Oj (s) functions as an ordinary derivative and obeys a chain rule; 

• a derivative with respect to an interaction coefficient J^i acting on a once or more than once primed quantity, 
such as m\ and ra", functions as an ordinary derivative and obeys the chain rule; 

• a derivative with respect to an interaction coefficient J^i acting on an unprimed quantity such as rrii must be 
treated apart, since this derivative will include a term taken on the tanh, which in turn will give a higher order 
correlation. 

These rules allow us to continue from what has already been computed and write 

dm i i /, 2xr e 

lind = (1 -m^dijdst 



3^ 772 " 

80^)30^)69^') ^ = ^-^H 3 ™*- 
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For the mixed terms we have similarly 



flT-lind = 0--m i )S ik m j 

UO j k 

d 2 m- 

dJ jk de l \s)d ^ md = ~ 2mi ^ ~ m i) s n s >,t8ikm' j + (1 - m 2 )S ik (l - (m^) 2 )<5ji*»,t-i 

- 2mj(l - m?)5i/<J S)t 5ife(l - (mj-) 2 )<Sj/j£ s / it _i 

- 2mj(l - m 2 )5ij/5 s ',t5ife(l - (m'j) 2 )6jid a> t-i 

+ (1 - m?)J ifc (-2mJ)(l - (m' j ) 2 )6 j i5 att - 1 5 jV 5 a , tt - 1 



and 



dJikdJi 



| ilH j = <J ifc (l - m?)(l - (m;-) 2 ),5 ifc m ; + (jfc) (Zm) 
-2m, (1 - m?)<5 ifc <y im (m' fc m^ + X ' fe J) 



dJ jk dJi m d8 n {s) 



lind = < M _2m i( 1 _ TO »)<*irA,t(l - K) )<5; fe m ; 

+^ fe (l - m 2 )(-2mj(l - (TOj) 2 )5j„(5 S!t _i(5( fc TO,' 
+<5 ife (l - m 2 )(l - (m' 3 ) 2 )S lk (l - (m';) 2 )S ln S s ^ 2 
+ (jk) <-> (Zm) 

+ 2(1 - m 2 )(3mi - l)^m^s,^ife^m < &k(t - l)<7 m (t - 1) > 
-2mj(l - m 2 )(5 ife <5 im (l - (m' fe ) 2 )(5fe„^ Sit _im^ 
-2m;(l - m 2 )5 lk 5 lm m' k (l - (m^) 2 )(5 m „(5 S!t _i 



-2wi(l - m^SikS. 



d6 n (s) 



where we use the correlation function Xfem =< v k {t)&m{t) > —m k m m . Its partial derivative with respect to an 
external field is always zero, and the last term in above therefore vanishes. The more cumbersome term is three 
derivatives with respect to interaction coefficients, which we can start from 



d 3 



OpgdJjkdJi 

r 



lind 



dJ, 



9 2 P( 



dJjkdJ, 
dP(a) 



tanh(-)+ 
(l-tanh 2 (-))<5 r 



+ (jk) o (Im) 
£ ff P(a)(-2tanh(-))(l-tanh 2 (.))4 



t ai(t - l)6ik<Tj(t - 1) 



(B7) 



Applying dJ pq gives (at least conceptually) eight terms. The term from acting on gj.^gj^ vanishes. The term from 
acting on tanh(-) in the first line gives a second derivative with respect to interaction coefficients of a magnetization. 
The terms from the second and the third line give combinations involving either second derivatives of a magnetization, 
or first derivatives of a correlation function. The terms from the last line are a third order correlation function and 
further first deritives of second order correlation functions. 
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Taking all together we can sum the contributions to 

Third order = -2(1 - m^)(3m^ - 1){A^ \t)f 



+ 2m,(l - m?)A l (1) (t)^(l - (mj) 2 ) J« A« (t 1) 
i 

(l-mj)Y,MAf\t)) 2 
i 

+ ^2(1 - m 2 )(3m 2 - l)^ 1} (i) £ JuJ im Xlm 

Im 

- 2m 2 (l - mf) M 1 ~ {m'iY)Af-\t - l)J m m' n 

In 

+ (1 - m 2 ) JmiJUl ~ «) 2 )(1 " (m';r)A^(t 2) + (m) O (0 

ml 

dXln{t-l)\ .(1) 



ln,js 



1 C1 2n t t f d Xln(t- 1)\ 

- o m Hl - m i ) J u J ni wj J Pq + circ. perm. 

A in,.,'* V / 

+ hi - m^)(3m^ - JuJniJ qi Xi„ q (B8) 



Inq 



where in the last line we have used xinq =< — m;)(cr„(i) — m n )(a q (t) — m g ) >. All the terms in above containing 

the first order terms A^ 1 ) vanish, the partial derivative terms of the second order correlation function with respect 
to external field vanish, and the last line is at least smaller than e 3 . The sole remaining terms hence come from 
the partial derivatives of second order correlation functions with respect to interaction parameters. These are model 
dependent, and are evaluated to non-zero for the sequential update rule in |18j . For the parallel update rule which we 
look at here they are however zero. The collection of terms (B8) therefore evaluates to zero. 



[1] D. Thouless, P. Anderson, and R. Palmer, Phil Mag 35, 593 (1977). 
[2] G. Parisi, Statistical field theory (Addison- Wesley, 1988). 
[3] T. Tanaka, Neural Computation 12, 1951 (2000). 

[4] J. S. Yedidia, W. T. Freeman, and Y. Weiss, Understanding Belief Propagation and its Generalizations (Science & Tech- 
nology Book, 2003), p. 239269. 

[5] M. Mezard and A. Montanari, Information, Physics and Computation (Oxford University Press, 2009). 

[6] J. Pearl, Probabilistic Reasoning in Intelligent Systems: Networks of Plausible Inference (Morgan Kaufmann, San Francisco, 
1988). 

[7] J. Pearl, Causality : Models, Reasoning, and Inference (Cambridge University Press, 2000). 
[8] M. Opper and O. Winther, Physical Review E 64, 056131 (2001). 

[9] F. R. Kschischang, B. J. Frey, and H.-A. Loeliger, IEEE Transactions on Information Theory 47, 498 (1998). 
[10] U. Frisch, Turbulence (Cambridge University Press, 1995), ISBN 0-521-45713-0. 
[11] N. van Kampen, Stochastic processes in physics and chemistry (Elsevier, 2007). 
[12] J. J. Hopfield, Proc. Natl. Acad. Sci. 79, 2554 (1982). 

[13] R. A. Blythe and A. J. McKane, Journal of Statistical Me chanics: Theory and Experiment 2007, P07018 (2007), URL 

|http : //stacks . iop . org/1742-5468/2007/i=07/a=P07018[ 
[14] B. Derrida, E. Gardner, and A. Zippelius, EPL (Europhysics Letters) 4, 167 (1987). 

[15] J. Hertz, G. Grinstein, and S. A. Solla, in Neural Networks for Computing AIP Conf Proc 151, edited by J. Denker (1986). 
[16] J. Hertz, G. Grinstein, and S. A. Solla, in Heidelberg Conference on Glassy Dynamics and Optimization, edited by I. Mor- 

genstern and J. L. van Hemmen (Springer Verlag, 1987). 
[17] A. Crisanti and H. Sompolinsky, Phys. Rev. A 37 (1988). 
[18] H. Kappen and J. J. Spanjers, Physical Review E 61, 5658 (2000). 

[19] Y. Kanoria and A. Montanari, Annals of Applied Probability (2011), in press, [arXiv:0907.0449[ . 

[20] I. Neri and D. Bolle, J. Stat. Mech. p. P08009 (2009). 

[21] E. Aurell and H. Mahmoudi, J. Stat. Mech. p. P04014 (2011). 

[22] E. Aurell and H. Mahmoudi, Communications in Theoretical Physics 56, 157 (2011). 



17 



[23] H. Sompolinsky and A. Zippelius, Phys. Rev. B 25, 6860 (1982). 

[24] A. Crisanti and H. Sompolinsky, Phys. Rev. A 36, 49224939 (1987). 

[25] A. C. C. Coolen, S. N. Laughton, and D. Sherrington, Phys. Rev. B 53, 8184 (1996). 

[26] H. J. Sommers, Phys. Rev. Lett. 58, 1268 (1987). 

[27] J. P. L. Hatchett, B. Wemmenhove, I. P. Castillo, T. Nikoletopoulos, N. S. Skantzos, and A. C. C. Coolen, J. Phys. A: 

Math. Gen. 37 (2004). 
[28] J. Hertz and Y. Roudi, Phys Rev Lett 106, 048702 (2011). 
[29] Y. Roudi and J. Hertz, J. Stat. Mech. p. P03031 (2011). 

[30] M. Mezard and J. Sakellariou, Journal of Statistical Mechanics: Theory and Experiment p. L07001 (2011). 
[31] D. Sherrington and S. Kirkpatrick, Phys. Rev. Lett 35, 1792 (1975). 
[32] A. Kholodenko, Journal of Statistical Physics 58, 355 (1990). 

[33] T. Plefka, Journal of Physics A: Mathematical and general 15, 1971 (1982); 

[34] J. Sakellariou, Y. Roudi, M. Mezard, and J. Hertz, Arxiv preprint [arXiv:1106. 04*52] (2011). 
[35] U. Marconi, A. Puglisi, L. Rondoni, and A. Vulpi ani, Physics Repo rts 461, 111 (2008). 
[36] J. Hertz, Y. Roudi, and J. Tyrcha, Arxiv preprint [arXiv:1106. 1752] (2011). 

[37] S. Amari, S. Ikeda, and H. Shimokawa, Information geometry and mean field approximation: the alpha-projection approach 

(MIT Press, 2001), pp. 241-257, ISBN 0-262-15054-9. 
[38] M. Mezard, G. Parisi, and M. Virasoro, Spin Glass Theory and Beyond (World Scientific, Singapore, 1987). 
[39] M. Mezard and G. Parisi, Eur. Phys. Journ. B 20 (2001). 



